Comparing the inspiral of irrotational and corotational binary neutron stars 
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We model the adiabatic inspiral of relativistic binary neutron stars in a quasi-equilibrium (QE) 
approximation, and compute the gravitational wavetrain from the late phase of the inspiral. We 
compare corotational and irrotational sequences and find a significant difference in the inspiral 
rate, which is almost entirely caused by differences in the binding energy. We also compare our 
results with those of a point-mass post-Newtonian calculation. We illustrate how the late inspiral 
wavetrain computed with our QE numerical scheme can be matched to the subsequent plunge and 
merger waveform calculated with a fully relativistic hydrodynamics code. 
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I. INTRODUCTION 

Binary neutron stars are among the most promising 
sources for the detection of gravitational waves by the 
laser interferometers currently under development, in- 
cluding LIGO, VIRGO, GEO and TAMA (see § for the 
first preliminary data analysis from the TAMA detector). 
Typical binary neutron stars are expected to complete 
about 16,000 orbital periods (corresponding to about 15 
minutes) while sweeping through LIGO's frequency band 
between 10 and 1000 Hz Accurate knowledge of the 
binary inspiral and the associated gravitational waveform 
is crucial to enhance our chances of detecting the signal 
in the noisy output of the laser interferometer. 

The evolution of binary neutron stars proceeds in dif- 
ferent stages. By far the longest is the initial quasi- 
equilibrium inspiral phase, during which the stars move 
in nearly circular orbits, while the separation between the 
stars decreases adiabatically as energy is carried away by 
gravitational radiation. The quasi-circular orbits become 
unstable at the innermost stable circular orbit (ISCO), 
where the inspiraling system enters a plunge and merger 
phase. The merger and coalescence of the stars happens 
on a dynamical timescale, and produces either a black 
hole or a larger neutron star, which may collapse to a 
black hole at a later time. The final stage of the evo- 
lution is the ringdown phase, during which the merged 
object settles down to equilibrium. 

The early inspiral phase, for large binary separations, 
can be modeled very accurately with post-Newtonian 
(PN) methods, e.g. H. Recently, the convergence of 
these PN expansion has been improved by the introduc- 
tion of Pade approximants M , and an alternative PN "ef- 
fective one-body" approach has also been suggested || 
(see also Q and references therein). It is generally ac- 
cepted that the plunge and merger phase has to be sim- 
ulated by means of a self-consistent, fully relativistic hy- 



drodynamics calculation (see Q for the first such sim- 
ulation). During the late ringdown phase, the merged 
object can be approximated as a distorted equilibrium 
object, so that perturbative techniques can be applied 
(e.g. §§)■ 

It is likely, however, that PN point-mass techniques 
break down somewhat outside of the ISCO, when finite- 
size and relativistic effects become important. It is hard 
to imagine that fully hydrodynamical numerical calcula- 
tions will be able to follow the inspiral reliably by begin- 
ning far beyond that point, through many orbital periods, 
all the way to the ISCO, followed by plunge and merger. 
Such calculations would accumulate significant amounts 
of numerical error and would be computationally pro- 
hibitive. This leaves a gap between the regimes that PN 
and fully numerical calculations can model. Filling this 
gap for the late inspiral, immediately prior to plunge, 
therefore requires an alternative approach (in the case 
of binary black holes, this problem has been called the 
"intermediate binary black hole" problem (]iof ). 

We have recently adopted a quasi-equilibrium (QE) 
approach to model the late inspiral This approach 
takes advantage of the fact that during the inspiral phase, 
the binaries evolve slowly on nearly circular orbits. The 
matter distribution can therefore be computed indepen- 
dently by assuming the stars to be in quasi-equilibrium 
(e.g. []l2 13 1 ) . These pre-determined matter profiles, 



which satisfy the hydrodynamic equations of quasi-static 
equilibrium, can then be inserted into Einstein's field 
equations as source terms and the hydrodynamic equa- 
tions do not have to be solved again (see the "hydro- 
without- hydro" approach discussed in [Q). For a given 
separation, we evolve the binary for about two orbital 
periods to determine the gravitational wave form and lu- 
minosity. We repeat the calculation for a discrete set of 
separations, which then allows us to construct the en- 
tire late inspiral together with its gravitational wave sig- 
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nal. The QE method is described in more detail in ]Tl| ] 
and in Section [I A below. We have also tested and cal- 



ibrated this method for a model problem in relativistic 
scalar gravitation [fil| , where the QE wavetrain could be 
compared with the full wavetrain calculated by an exact 
numerical integration. 

In we presented a prototype calculation for a coro- 
tational binary sequence. It is unlikely, however, that 
binary neutron stars can maintain corotation during in- 
spiral Jl6[ . In this paper, we therefore compare with the 
more realistic case of irrotational binary neutron stars. 
By comparing the results we thereby explore the sensi- 
tivity of the wavetrain to the internal fluid velocities. We 
discuss the criterion which determines the point of break- 
down of the QE approximation just outside the ISCO. We 
thus establish the innermost separation at which a QE 
calculation can provide initial data for future dynamical 
simulations of the plunge and merger. Finally, we com- 
pare our computational results for the late inspiral with 
a point-mass PN calculation of the same epoch. 

The paper is organized as follows. In Section || we 
describe our numerical methods and resulting waveforms 
(Section [II A ), present a first order PN formalism (Sec- 
tion [II B| ) , explain the construction of an inspiral wave- 
train (Section [I C), and finally compare corotational and 
irrotational sequences, in both the numerical and PN ap- 



proach (Section II D). In Section III we illustrate how 



an entire wavetrain, spanning both the late inspiral and 
the plunge and merger phases, can be constructed by 
matching QE and fully dynamical results. We briefly 



summarize our results in Section IV 



II. 



IRROTATIONAL AND COROTATIONAL 
BINARY INSPIRAL 

A. Numerical results 



We adopt a quasi-equilibrium (QE) approximation to 
study the adiabatic, late inspiral of binary neutron stars, 
taking advantage of the fact that in this regime the or- 
bital decay timescale is much longer than the orbital pe- 
riod, and, in addition, the gravitational radiation has 
negligible effect on the structure of the stars (see for 
details). In this approximation, the inspiral is modeled 
as a sequence of "snapshots" of binaries in quasi-circular 
orbit. In this paper we generalize the results of for 
corotational binaries (based on the models of |12|l ) to ir- 
rotational binaries (based on the models of fl3fl) and an- 
alyze the effect of spin on the inspiral. The snapshots 
used were constructed under the assumptions that the 
metric is asymptotically and conformally flat, although 
the metric ceases to be conformally flat as we evolve it 
to generate waveforms. These models provide the data 
to construct binding energy curves, which give the to- 



tal mass-energy M(r) as a function of separation r for 
circular binaries of fixed rest (baryonic) mass. 

We use these binary models to provide the matter 
source terms in the Einstein field equations, which we 
integrate numerically with the formalism and code de- 
scribed in Q (see also |Q) and used in 0|. For bi- 
nary separations where the QE approximation is valid, 
the hydrodynamic equations do not need to be resolved 
after the initial time . Instead, we freeze the matter vari- 
ables at their equilibrium values and move the stars on 
circular orbits while allowing the metric to evolve dy- 
namically. We typically evolve the binary for about two 
orbital periods for each orbital separation, which is suf- 
ficient to read off the gravitational waveform h and the 
gravitational wave luminosity —dM/dt. The initial grav- 
itational wave signal is corrupted by noise contained in 
the approximate initial data, which quickly propagates 
off the numerical grid leaving a periodic wave signal from 
the binary. Repeating the integration for various binary 
separations r yields h(r) and dM/dt(r) as a function of 
separation. The inspiral rate can then be found by com- 
bining the luminosity with the derivative of the binding 
energy with respect to separation, 



dr _ dM/dt 
di ~ dM/dr' 



(1) 



Integrating this equation yields the separation r as a 
function of time t. The gravitational wavetrain h(t) 
can finally be constructed by using r(t) to fit the dis- 
crete waveforms h(r) together in a smooth and contin- 
uous fashion. This QE approach has also been tested 
and calibrated for a model problem in relativistic scalar 
gravitation ( Jl5|). 

Both the corotational and the irrotational sequences 
terminate at some innermost circular orbit (ICO), at 
which the quasi-equilibrium approximation breaks down. 
For the corotational sequences, this ICO corresponds to 
the innermost stable circular orbit (ISCO) at which the 
equilibrium orbits become unstable, and at which the 
plunge and merger of the binary set in. Typical irrota- 
tional equilibrium sequences terminate before they reach 
the ISCO, when they form cusp-like structures at the in- 
ner edges of the stars Jl3],[l9]]. The cusps indicate the 
onset of mass overflow. The termination of the equilib- 
rium sequence indicates the break-down of the QE ap- 
proximation, since no equilibrium object exists with this 
separation. For both the corotational and the irrotational 
sequences, we denote the mass-energy at the innermost 
circular orbit as Mico- 

By adopting QE matter sources in our evolution cal- 
culations, we assume that these remain stationary in a 
co-rotating frame of reference over an orbital period. As a 
necessary condition for this to be true, we have to ensure 
that the binary does not reach the ICO in approximately 
one orbital period P. Since the energy emitted by the 
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TABLE I. Irrotational Sequence 



ZA a 


M/M b 


fiorb-Afo 


AM C 


dM/dt 


r/M d 


dr/dt 


iVc y c e 




0.500 


0.95261 


0.00587 


0.0190 


4.94e-09 


14.7 


2.21e-05 


0.0 


0.001582 


0.474 


0.95237 


0.00648 


0.0219 


8.022e-09 


13.7 


2.79e-05 


77.1 


0.00250 


0.444 


0.95210 


0.00722 


0.0234 


1.14e-08 


12.7 


3.73e-05 


144.9 


0.00351 


0.412 


0.95177 


0.00808 


0.0261 


1.78e-08 


11.8 


5.12e-05 


198.3 


0.00552 


0.375 


0.95141 


0.00913 


0.0290 


2.80e-08 


10.8 


7.37e-05 


242.4 


0.00898 


0.000 


0.95100 


0.0104 


0.0335 


A O A ~ no 

4.o4e-U8 


n on 
9.89 


n nnm 1 n 
0.000110 


077 n 
I 1 1 .U 


0.0168 


0.286 


0.95052 


0.0119 


0.0379 


8.20e-08 


8.99 


0.000169 


302.8 


0.0343 


0.231 


0.95001 


0.0138 


0.0440 


1.47e-07 


8.10 


0.000277 


322.1 


0.0896 


0.167 


0.94949 


0.0159 


0.0525 


2.79e-07 


7.28 


0.000563 


334.9 


0.483 


* 0.130 


0.94929 


0.0169 


0.0569 


3.69e-07 


6.91 


0.00135 


338.5 


4.90 


* 0.111 


0.94926 


0.0170 


0.0580 


3.87e-07 


6.78 


0.0121 


339.1 


oo 


a za is 


defined to be the ratio of the separation between the innermost points and the outermost points on 


the stars. 


b M is 


the total rest mass of the system, i.e., twice the rest mass of an individual star. 






c A is the amplitude of 


gravitational waves on 


the rotation axis 










d r is denned as the average of the coordinate distance from the 


origin of the nearest and farthest points on a neutron 


star. Thus, it is the coordinate radius of the "center" in n = 1 units and the metric and coordinate system 


of mm 


e -^cyc is the number of 


gravitational wave cycles from our initial configuration (i.e. that of largest r) to the indicated 


separation. 
















* QE approximation is 


invalid: £ > 1. 




















TABLE II. Corotational Sequence 








ZA 


M/M 


S7 or b Mo 


AMo 


dM/dt 


r/Mo 


dr/dt 


Acyc 




0.475 


0.952098 


0.00656 


0.0217 


8.06e-09 


13.7 


4.57e-05 


0.0 


0.00343 


0.45 


0.95190 


0.00715 


0.0236 


1.13e-08 


12.9 


3.66e-05 


46.6 


0.00486 


0.400 


0.951493 


0.00843 


0.0270 


2.07e-08 


11.4 


6.40e-05 


128.1 


0.00938 


0.375 


0.95131 


0.00922 


0.0299 


3.04e-08 


10.8 


9.32e-05 


152.0 


0.01416 


0.350 


0.951106 


0.0102 


0.0341 


4.82e-08 


10.1 


0.000130 


168.8 


0.0236 


0.300 


0.950752 


0.0114 


0.0368 


7.06e-08 


9.27 


0.000242 


189.8 


0.0429 


0.275 


0.950569 


0.0122 


0.0393 


9.20e-08 


8.84 


0.000320 


196.5 


0.0656 


0.250 


0.950471 


0.0130 


0.0415 


1.16e-07 


8.45 


0.000422 


201.6 


0.09014 


0.225 


0.950329 


0.0138 


0.0446 


1.52e-07 


8.09 


0.000557 


205.5 


0.143 


0.200 


0.950228 


0.0146 


0.0477 


1.93e-07 


7.78 


0.000744 


208.4 


0.219 


0.175 


0.950093 


0.0154 


0.0500 


2.36e-07 


7.47 


0.00101 


210.6 


0.393 


0.15 


0.950031 


0.0161 


0.0540 


3.02e-07 


7.21 


0.00140 


212.2 


0.645 


0.125 


0.950024 


0.0168 


0.0565 


3.59e-07 


6.98 


0.00204 


213.3 


0.761 


* 0.100 


0.949935 


0.0174 


0.0595 


4.31e-07 


6.75 


0.00312 


214.1 


1.77 



binary per period is (dM/dt)P, this condition holds as 
long as the dimcnsionlcss ratio 



\dM/dt\P 
= M - Mico 



(2) 



is less than unity. Identifying the orbital timescale t rb 
with P and the gravitational wave inspiral timescale tew 
with M / \dM / dt\ shows that the condition £ < 1 is iden- 
tical to 



M - Mico t mh 



M 



tew 



(3) 



Eq. (JsJ) states that, for infinitesimally small torbAcWj we 
can apply the QE approximation arbitrarily close to the 
ICO. However, if £ > 1, the binaries will overshoot the 
innermost equilibrium orbit from their current separation 



in less than one orbital period, causing a breakdown in 
quasi-equilibrium. 

The above criterion is not only useful for determining 
the limits of the QE approximation, but also for choosing 
the separation for assigning initial data for dynamical 
simulations of the plunge and merger. In order to reduce 
computational resources and accumulation of numerical 
error, one would like to impose such initial data as close 
to the ICO as possible. Very close to the ICO, however, 
the approximations of the QE approach break down, as 
indicated by Eq. (0). The separation at which £ sa 1 
is therefore a suitable compromise for terminating the 
QE approach, and imposing initial data for dynamical 
simulations. 

In Tables I and II we summarize our numerical re- 
sults for irrotational and corotational binaries. For both 
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B. Post-Newtonian treatment 
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t/M 

FIG. 1. The L2 norm of the Hamiltonian constraint and of 
the x and y components of the momentum constraint for coro- 
tating binaries with separation parameter za =0.1 orbiting 
about the z-axis. The Hamiltonian constraint is normalized 
to the initial value of 16nM, the momentum constraint to the 
L2 norm of the initial value of 8irS x , where S x is the domi- 
nant component of the initial matter current. The L2 norm 
is defined as L2(f) = y/Yj! 



cases the stars are modeled as polytropes with pressure 
P = K/9Q +1 ^ n , where po is the rest-mass density, n is the 
polytropic index, and k is a constant. We set n = 1 and 
non-dimensionalize our results by setting k = G = c = 1 . 
In these units, the individual stars have a rest mass 
too = 0.1, which is about 55% of the maximum allowed 
rest mass of isolated, spherical nonrotating stars, and 
corresponds to a compaction of rrioo/Roc — 0.088 in iso- 
lation. Here and are the mass-energy and areal 
radius (i.e. radius defined by the area of a constant-r 
shell) of the individual (spherical) stars at infinite binary 
separation. Our results apply to arbitrary mass stars 
with this compaction and equation of state. For both 
cases, we study models up to and past £ = 1. In Figure 
1, we verify that our evolution code satisfies the Hamil- 
tonian and momentum constraints. This is a nontrivial 
test of our method: The matter profiles are determined 
by solving the constraint equations in the QE approxi- 
mation. The fact that they remain satisfied during the 
field evolution demonstrates that the spacetime remains 
close to quasi-equilibrium, even as gravitational radiation 
is generated and the metric is determined dynamically. 



We adopt the formalism developed in p(I to study the 
binary inspiral in a first order PN approximation (note 
that PN techniques have been extended to at least third 
order; see, e.g. Q). In particular, we consider a system of 
two point masses in circular orbit with binary separation 
r, total orbital angular momentum Ln, and total spin 
angular momentum S. In isolation, each point mass has 
a mass-energy TOqo. From p0[ | we then find the orbital 
frequency 



- f* h 11 ^ 5 t ^ 



the total mass-energy 
and the luminosity 



(4) 



(5) 



= -1 (_) 5 (l - + 47r(-) 3/2 

dt 5 l r j 1 42 r K r J 



4(-) 3/2 L N 

3 r' 



s). 



Here we identify 



2TOoc 

2(to c 



±m 2 1 



irrotational 
corotational 



(6) 



(7) 



for the total mass-energy of each star in isolation. To 



compare with the binaries of Section II A , we choose 



for the total orbital angular momentum and 



S = 







irrotational 



2/£! or be z corotational 



(8) 



(9) 



for the spin angular momentum. 

In the above, r is the separation of the point masses in 
harmonic coordinates, and the moment of inertia of each 
mass, /, is calculated for a slowly rotating (nearly TOV) 
relativistic star Note that, for the corotational case, 
H, and O or b must be solved simultaneously to PN order. 



C. Constructing Wavetrains 

We extract gravitational wave data from our QE simu- 
lations by matching the numerical data near the bound- 
aries of our grids to the I = 2, m = ±2 Moncrief |22j ] 
variables R22+ and R22- (see Appendix A for details). 
From these, the energy luminosity is given by 



dM 
~dt 



32tt 



[(d t R22+) 2 + (d t R22-) 2 ]- 



(10) 
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FIG. 2. The total mass-energy M as a function of orbital 
frequency for corotational and irrotational models. The cross 
marks the ISO, the asterisk the location at which £ = 1. 



To generate the wavetrain, we also need to determine 
the wave amplitude. Both the QE and PN models pro- 
vide more information than just the quadrupole ampli- 
tudes, but we have found [jllj that the amplitudes of other 
modes are much smaller than the I = 2, m = ±2 modes, 
and we ignore them below. The asymptotic amplitudes 
of h+ and h x for these modes are identical along the axis 
of rotation, and are called A below. This amplitude is 
related to the luminosity by the quadrupole relation 



A 



1 



GW 



dM 



10- 



dt ' 



(11) 



where f^cw = 2f2 or b is the frequency of the waves. This 
gives A(r), which can be converted to A(t) by solving 
Eq. (0). Then the complete wavetrain as seen by an 
observer on the rotation axis at a distance r s from the 
source is given by 



r s h(t) = A(t) cos( / tt GW (t')dt') 



(12) 



D. Comparing irrotational and corotational inspiral 

In Figure 2, we show the total mass-energy M of irro- 
tational and corotational inspiral sequences (along which 
the total rest mass M is conserved). We parameter- 
ize the binary separation by fi or b, the orbital angular 
frequency as measured by a distant observer, since it is 
an invariant. Note that infinite separation corresponds 
to r2 or b = 0. The PN and QE results agree quite well 



— 



4x10- 



3x10- 



2xl0- 7 



10" 





0.005 



_ 1 1 1 


i I 1 I I i I 




" PN 


; 

irrot. / 




- 

- PN 


/ 

corot. A 


- 
— 


- 


fm 


- 


■ QE 


irrot J 








a QE 


corot / A 






/ ft 






/ A 






/ A 






/ ■ 






/ A 












* 






1 1 1 1 1 1 


1 



0.0 



0.015 



"orb M 



FIG. 3. The gravitational wave luminosity dM/dt as a 
function of orbital frequency for corotational and irrotational 
models. The solid and dotted curves nearly coincide. 
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FIG. 4. Gravitational wave amplitude along the axis of ro- 
tation in the QE, PN and Newtonian quadrupole approxima- 
tion. The PN curves for the two cases, which are based on 
Eq. (Illl), are hardly distinguishable. 
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FIG. 5. The final two hundred cycles of the inspiral wave- 
form for corotating (top) and irrotating (bottom) binaries. 



for large separations, while for small separations we find 
growing deviations between the two approaches. This is 
not surprising, since at small separations both tidal in- 
teractions and relativistic effects play an increasingly im- 
portant role. Finite size effects are not taken into account 
by the point-mass calculation of the PN approach. The 
inclusion of 2 and 3PN terms will not improve the agree- 
ment substantially, due to the dominance of finite size 
effects ^,^4|. For small separations, the mass-energy, 
and hence the binding energy, of irrotational and coro- 
tational sequences is quite different. This effect is due 
to the extra spin rotational energy of the corotating bi- 
naries, which causes their total energy to be larger than 
that of the irrotational binaries. As we will see below, 
this change in M, and hence dM/dr, dominates the dif- 
ference in the inspiral rate dr/dt. 

In Figure 3, we plot dM/dt for each model. Note that 
the gravitational wave luminosity dM/dt is very simi- 
lar for the irrotational and corotational sequence in ei- 
ther the QE or the PN approximation. In Figure 4, 
we show the gravitational wave amplitudes as measured 
along the axis of rotation and compare two different nu- 
merical results with the PN values. From our numerical 



binary models, we predict the gravitational wave ampli- 
tude, both from the Newtonian quadrupole approxima- 
tion, where the quadrupole moments are computed by 
numerical quadrature of the density and velocity fields, 
and by dynamically evolving the (relativistic) gravita- 
tional fields and reading off the wave amplitudes on the 
outer edge of the numerical grid, as described in section 
IIC. 

For large separations, the quadrupole and PN results 
agree quite well, while the QE results differ by ~ 30%. 
We attribute this difference to the inaccurate boundary 
treatment in our dynamical evolution (see ) . Gravita- 
tional radiation can only be read off accurately in the far 
(wave) zone, at separations greater than a gravitational 
wavelength Aqw- I n our simulations, however, we are 
restricted to numerical grids which only extend to a frac- 
tion of a gravitational wavelength, r max ~ 0.1 — 0.35Aqw- 
For large binary separations, the frequency SIgw is small, 
and hence Aqw is large. For fixed outer boundary 
r max this implies that r max /Xow is smaller for larger 
binary separations, causing a more inaccurate wave ex- 
traction. As the separation becomes smaller, the QE 
and quadrupole results agree increasingly well. Simulta- 
neously we expect the PN point-mass treatment to be- 
come increasingly inaccurate for small separations, since 
it neglects finite size effects. This trend is reflected in 
the growing deviations between the PN and quadrupole 
results in Figure 4. 

The differences between the irrotational and corota- 
tional wave amplitudes are very small in all three ap- 
proaches. This is quite intuitive, since the gravitational 
wave emission is dominated by the matter density dis- 
tribution, which is fairly similar for the two sequences, 
while matter current distributions play a less important 
role. The poor handling of the outer boundaries in the 
QE approximation affects both sequences in the same 
systematic way, which allows us to make a meaningful 
comparison. 

The findings of this comparison with regard to the in- 
spiral rate can be anticipated from Eq. p). The numera- 
tor on the right hand side, dM/dt is very similar for irro- 
tational and corotational sequences, while the denomina- 
tor, dM/dr is different. The latter therefore causes the 
inspiral rate | dr/dt | to be smaller for irrotational than 
for corotational binaries. This effect can be understood 
very easily even by a crude Newtonian argument, where 
the total mass-energy of the star can be written as a sum 



U 2r 



-2-m 2 

2 



M - 



2r 



(13) 



In Eq. (|T^), the first term on the right-hand side is the 
energy of a spherical polytrope (rest + internal + grav- 
itational potential energy), the second term is the com- 
bined orbital kinetic and gravitational potential energy 



G 



for circular equilibrium, and the third term, in square 
brackets, denotes the spin kinetic energy. Note that the 
term in the square brackets applies only to the corota- 
tional sequence. For a given separation, the energy of a 
corotational binary is thus larger than that of an irrota- 
tional binary by the spin energy of the individual stars. 
Taking a derivative with respect to r yields 



dM corot dM mot M 
dr dr r 4 



(14) 



Inserting this into Eq. (|l]), and using the result that 
dM/dt is nearly identical for both sequences, immedi- 
ately shows that 



> 



(15) 



irrot 



Correspondingly, the frequency of the gravitational wave 
train increases at a faster rate for corotational binaries, 
which we illustrate in Figure 5. 



III. 



COMPUTING THE FULL WAVETRAIN 
THROUGH COALESCENCE 



As we have argued in Section [I A, the QE approach 
could provide accurate numerical initial data for fully 
relativistic hydrodynamics simulations of the plunge and 
merger at a separation for which £ <; 1. Starting at 
this separation guarantees that the binary reaches the dy- 
namical instability at the ICO within about one orbital 
period. This approach provide intrinsically more accu- 
rate initial data than those imposed at the ICO (where 
the QE approximation breaks down). Moreover, the QE 
approach allows noise contained in the initial data to be 
radiated away, thereby producing more accurate equilib- 
rium data. In the QE approach, the initial noise can be 
easily distinguished from the later periodic binary sig- 
nal, whereas in a dynamical plunge simulation it may 
be more difficult to distinguish the propagation of initial 
noise from the intrinsically short-lived coalescence wave- 
form. 

In Figure 6, we provide an illustration of what the com- 
bined late inspiral-plunge wavetrain might look like. We 
compute the late inspiral wavetrain using the QE method 
for an irrotational binary sequence with stars of rest mass 
too = 0.146, which corresponds to 81% of the maximum 
allowed rest mass of an isolated, nonrotating star, and 
to a compaction of rrioo/Roo =0.14 in isolation p5| . Wc 
then match this to the full hydro dynamical waveform for 
the plunge and merger of this same binary as computed 
by Shibata and Uryu j7lp7J ■ In order to get a continuous 
waveform, we adopt the same wave extraction rules as in 
0, namely, h x is read off the rotation axis at the edge of 
the numerical grid, located at r — 0.35Aqw- The contin- 
uous wavetrain is constructed in the following manner: 
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FIG. 6. A match of the total late QE inspiral wavetrain to 
the plunge and merger waveform, both constructed from nu- 
merical simulations in general relativity. Here we connect the 
last 13 orbits of the inspiral, computed in the QE approxima- 
tion, to the plunge calculation of Q. The ICO is reached at 
t/Mo =1330, around which point we pass from the adiabatic 
inspiral to the plunge and merger phase. Here h x = x^7xy, 
where r 3 is the distance from the source and 7y is the confor- 
mal 3-metric on the rotation axis. This coalescence scenario 
leads to merger and immediate black hole formation. For 
Mo = 2 x 1.5Mq and r 3 — 100 Mpc, the maximum amplitude 
of the metric perturbation is ~ 5 x 10~ 21 . 



From t = to t = 1321M at the ICO, we use QE data. 
A smooth (sine) function with the appropriate gravita- 
tional wave frequency and amplitude is used to match 
the QE data to the plunge data. The raw plunge data 
given m begins at t = 1367M . 

The location and nature of the transition from quasi- 
periodic inspiral to rapid plunge is of great astrophysical 
interest because it is sensitive to details of the binary neu- 
tron star system, in particular the equation of state and 
rotation rate p^j. As illustrated above, the QE method, 
in conjunction with a relativistic hydrodynamics code, 
could allow one to obtain the complete wavetrain from 
late inspiral through plunge and merger, and at a very 
reasonable computational cost. 



IV. SUMMARY 

We study the late inspiral of binary neutron stars and, 
by comparing corotational and irrotational sequences, ex- 
plore the effect of internal fluid motions on the inspiral 
rate. We find the gravitational wave luminosity is very 
similar for any given radius. However, the inspiral of 
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corotational binaries proceeds faster. 

We also compare our numerical QE results with a first- 
order PN calculation. For large binary separations, the 
poor outer boundary conditions in the numerical calcu- 
lations introduce significant deviations, while for small 
separations finite-size effects, which are not taken into 
account in the PN calculations, can no longer be ne- 
glected. Future improvements in wave extraction and 
outer boundary conditions (see, e.g. p9[) should reduce 
the discrepancy. These improvements remain a long-term 
goal for numerical simulations of binary inspiral. 

Lastly, we illustrate how QE and fully dynamical re- 
sults could be matched to produce the entire gravita- 
tional wavetrain, covering the late inspiral through the 
plunge and merger phases. We emphasize that this pa- 
per presents only a prototype calculation showing how, 
in principle, the continuous wavetrain can be assembled 
from numerical simulations. 



(Y22 +y 2 -2)/v / 2, and R22-, which is associated with the 
mode having angular dependence (Y22 — Y 2 _ 2 )/i\/2- 

As argued in ]3(J], our gauge conditions approach the 
TT gauge as r becomes large. In this gauge, the two 
gravitational wave polarizations can be constructed from 



lm 



h x = r^G hn X lm / s'm( 

I in 



(A2) 



The Moncrief functions only depend on r and t, while the 
angular dependence is contained in the Wi m and Xi m . In 
the TT gauge, all Moncrief variables vanish except for the 
Gim, so that the linear combination for the Rf m reduces 
to 



2±2 



2V12G 2±2 



(A3) 



for the dominant modes I — 2, m — ±2. Inserting these 
into eq. (|A3) we then construct 
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APPENDIX A: GRAVITATIONAL WAVE 
EXTRACTION 

The two gravitational wave polarizations can be ex- 
pressed as 



h+ = "-J r|-(l + cos 2 9)(cos(2cj))R 22+ + sin(20)i? 22 _) 
2 V 1D7T 



/' = cos<9(- sin(2</>)i? 22+ + cos(2</>)i? 22 _) 

(A4) 

up to an arbitrary phase. Finally, the amplitude A on 
the z-axis (i.e. 9 = 0) is given by 



A = r X /—^(R 22+ y + (R 22 -y. 



(A5) 



As a check of the above, one can, by making explicit 
the time dependence in (A4) (4> — > 4> — O, or ^t), use 
Eq. ( |Io| ) to derive the luminosity-amplitude relationship 
for quadrupole (e.g. point mass) systems, Eq. (|TT|) . The 
amplitude of the waveform h(t) in Eq. ( |l2| ) is then equal 
to A(t)/r. (The 1/r falloff is restored.) 



TT 



u TT 



Irhrh > 



(Al) 



where the superscript TT denotes transverse-traceless 
projection, and 7jj is the spatial metric in orthonormal 
spherical-polar coordinates. Note that the amplitudes 
of h + and h y are independent of r at large r. We 
determine h + and h x from the Moncrief |p2|| variables 
Gi„i, hum, H 2 i m and Ki m , which we compute by quadra- 
tures of our metric over a spherical shell near the edge 
of our grid (see f|(J). From these we find the gauge- 
invariant linear combinations Rf m . In our code, we work 
with real, orthonormal combinations of spherical har- 
monics, and so we compute the amplitudes R 22 +, which 
is associated with the mode having angular dependence 
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